##################################################################

# Gautam Nair
# gautam.nair@yale.edu
# Misperceptions of Relative Affluence and Support for International Transfers
# Make Figure 2: Average Treatment Effects

##################################################################

# setwd("")

##################################################################
# loading packages
##################################################################
library(Hmisc)
library(foreign)
library(sandwich)
library(lmtest)
library(numDeriv)
library(stargazer)
library(ggplot2)
library(plyr)
library(gridExtra)
library(ri)
library(dplyr)
library(plyr)
library(scales)

##################################################################

rm(list=ls())
data.working <- readRDS("d_r_cleaned_analysis_dataset.Rda")

##################################################################
# Treatment Effects Figures
##################################################################

tbl_df(data.working)

outcome.vars <- c(
"q13ab.3",
"q13ab.3.positive",
"q8.toolittle",
"q11.yes",
"q8.toolittle.1.q11.yes.1",
"q10.reversed"
)

title.list <- c(
"Mean Percentage of $20 Payment \nDonated to International Charity",
"Proportion of Respondents Making a \n Positive Donation to International Charity",
"Proportion of Respondents Supporting\nIncreases in Foreign Economic Assistance",
"Proportion Supporting Lower \n Trade Protections for Domestic Agriculture",
"Prop. Supporting Increased Foreign \nAssistance and Lower Trade Protections",
"Mean Support for \n Domestic Redistribution (1-7)"
)

file.list <- c(
"int_charity_amount",
"int_charity_prop",
"foreign_aid_too_little",
"agri_cut_protections",
"aid_too_little_agri_cut",
"domestic_redistribution"
)

for(i in 1:length(outcome.vars)){
	data.working$Y <- data.working[, outcome.vars[i]]
	data.working$D <- data.working$group
	outcome <- data.working
	outcome <- group_by(outcome, D)
	# different versions based on whether outcome is dummy or continuous
	# dummy
	if(max(data.working$Y, na.rm=TRUE)==1){
		outcome <- summarize(outcome, 
		count=n(),
		y_bar=mean(Y, na.rm=TRUE),
		y_bar_100=format(round((mean(Y, na.rm=TRUE)*100), digits=1),nsmall=1),
		se=sqrt(y_bar*(1-y_bar)/n()),
		upper95 = y_bar + 1.96*se,
		lower95= y_bar - 1.96*se,
		upper90 = y_bar + 1.65*se,
		lower90= y_bar - 1.65*se)
	}
	# continuous
	if(max(data.working$Y, na.rm=TRUE)>1){
		outcome <- summarize(outcome, 
		count=n(),
        y_bar=mean(Y, na.rm=TRUE),
        y_bar_100=format(round((mean(Y, na.rm=TRUE)), digits=1),nsmall=1),
        se=sd(Y, na.rm=TRUE)/sqrt(n()),
        upper95 = y_bar + 1.96*se,
        lower95= y_bar - 1.96*se,
        upper90 = y_bar + 1.65*se,
        lower90= y_bar - 1.65*se)
	}	
	outcome[1,1] <- 2
	outcome[2,1] <- 4
	outcome[3,1] <- 6
	ate.plot <- ggplot(outcome, aes(x=D, y=y_bar, label=y_bar_100)) + geom_text(size=2.5, hjust=-0.8) +
						geom_pointrange(aes(ymin=lower90, ymax=upper90), size=0.8, color="navyblue") + 
						geom_pointrange(aes(ymin=lower95, ymax=upper95), size=0.4, color="navyblue") + 
						ggtitle(title.list[i]) + 
						xlab("") + ylab("")  +
						theme_bw(base_size = 9) + 
						theme(plot.title = element_text(size = 9), panel.grid.minor= element_blank(), axis.text=element_text(size=9), 
						axis.title=element_text(size=9)) +
      					scale_x_continuous(name="", breaks=c(2,4,6), labels=c("Control","Salience","Information"), limits=c(1,7)) +
						theme(plot.margin = unit(c(0.1,0.1,0.1,0.1), "cm"))
	# dummy					
	if(max(data.working$Y, na.rm=TRUE)==1) {
		ate.plot <- ate.plot + scale_y_continuous(labels = percent_format())
	}
	tempname <- paste("tf", "f", "02", i,  file.list[i], ".png",  sep="_")
	ggsave(ate.plot, filename=tempname, width=3,height=3)
}						
						
